Upper limit on the flux of photons with energies above 10 19 eV using 
Telescope Array surface detector 
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We search for ultra-high energy photons by analyzing geometrical properties of shower fronts of 
events registered by the Telescope Array surface detector. By making use of an event-by-event 
statistical method, we derive upper limits on the absolute flux of primary photons with energies 
above 10 19 , 10 19 ' 5 and 10 20 eV based on the first three years of data taken. 
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PACS numbers: 98.70.Sa, 96.50.sb, 96.50.sd 



I. INTRODUCTION 



hadron-induced EAS 



/ 



gamma-induced EAS 



The Telescope Array (TA) experiment pQ is a hybrid 
ultra-high energy (UHE) cosmic ray detector covering 
about 700 km 2 in central Utah, USA. It is composed 
of a Surface Detector (SD) array and three Fluorescence 
Detector (FD) stations. The TA SD array consists of 507 
plastic scintillator detectors on a square grid with 1.2 km 
spacing [3J. They each contain two layers of 1.2 cm thick 
plastic scintillator 3 m 2 in area. The three FD stations [3] 
contain a total of 38 telescopes overlooking the air space 
above the array of scintillator detectors. The purpose of 
this paper is to present the photon search capabilities of 
the Telescope Array surface detector and to search for 
primary photons in the cosmic ray flux. We place the 
limits on the integral flux of photons for energies greater 
than E , where E takes values 10 19 , 10 19 5 and 10 20 eV. 

At present there is no experimental evidence for 
primary UHE photons. However, several limits on 
the photon flux have been set by independent experi- 
ments. These include Haverah Park g], AG AS A [5], 
Yakutsk [6j [7] (see also reanalyses of the AGASA [8] 
and AGASA+Yakutsk [9] data at energies greater than 
10 20 eV) and the Pierre Auger Observatory jTUHH]. 

Photon limits may be used to constrain the parameters 
of top-down models [13]. In future the photon searches 
may be used to assess parameters of astrophysical sources 
in the Greisen-Zatsepin-Kuzmin Q31 [T5] cut-off scenario 
which predicts photons as ever present secondaries. If 
UHE photons are observed, they will be a supporting 
evidence for GZK nature of the spectrum break at the 
highest energies observed by HiRes [IB], Pierre Auger 
Observatory [T7] and TA [TS]. Photon flux is sensitive 
to the mass composition of cosmic rays and hence may 
be used as a probe of the latter (THl HO] . The results of 
the photon search also constrain parameters of Lorentz 
invariance violation [2THSB"] . Finally, photons with ener- 
gies greater than ~ 10 18 eV could be responsible for CR 
events correlated with BL Lac type objects on an angu- 
lar scale significantly smaller than the expected deflection 
for protons in cosmic magnetic fields, suggesting neutral 
primaries [26, 27] (see Ref. [28] for a possible mechanism). 

Since the TA detectors are composed of thin scintil- 
lators, they respond equally to the muon and electro- 
magnetic components of the extensive air shower and 
are therefore sensitive to showers induced by photon pri- 
maries (see e.g. Ref. [5§] for discussion). We use the 
shower front curvature as a Composition-sensitive pa- 
rameter (C-observable) and a modification of an event- 
by-event statistical method |30j to constrain the photon 
integral flux above the given energy. For the Energy- 
sensitive parameter (E-observable) , we use the scintilla- 
tor signal density at 800 m core distance S = Sgrjo- The 
comparison of event-by-event statistical method with the 
"photon median" method [TT] is presented. 
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FIG. 1: Illustrative view of hadron- and gamma-induced 
showers. Gamma-induced shower is deeper due to smaller 
cross-section of the first interaction. Moreover, the hadronic 
cascade is secondary with respect to electromagnetic in 
photon-induced showers. The latter contains fewer muons 
(shown in red) and have larger curvature of the shower front. 



II. SIMULATIONS 

Extensive Air Showers (EAS) induced by photon pri- 
maries differ significantly from hadron-induced events 
(see e.g. [3T] for a review). Photon induced showers con- 
tain fewer muons and have a deeper shower maximum 
when compared to hadronic showers. The latter results 
in the shower front having more curvature at the sur- 
face as illustrated in Figure [T] At the highest energies, 
there are two competing effects responsible for the diver- 
sity of showers induced by photon primaries. First, the 
electromagnetic cross-section is suppressed at energies, 
E > 10 19 eV due to the Landau, Pomeranchuk [55] and 
Migdal [33] (LPM) effect. The LPM effect delays the first 
interaction so that the shower arrives at ground level un- 
derdeveloped. The second effect is e ± pair production 
which is due to photon interaction with the geomagnetic 
field above the atmosphere. Secondary electrons produce 
gamma rays by synchrotron radiation generating a cas- 
cade in the geomagnetic field. The probability of photon 
conversion is a function of photon energy and the perpen- 
dicular component of geomagnetic field [31] . The shower 
development therefore depends on both the zenith and 
azimuthal angles of the photon arrival direction. 

The event-by-event method [30 requires a set of simu- 
lated photon- induced showers for the analysis of each real 
shower. We simulate the library of these showers with 
different primary energies and arrival directions. For the 
highest energy candidates (events which may be induced 
by a photon with primary energy greater than 10 19 ' 5 eV) 
we simulate individual sets of showers with fixed zenith 
and azimuthal angles. At these energies the shower de- 
velopment becomes azimuth-angle dependent due to the 
photon cascading in the geomagnetic field |31j . 

We use CORSIKA [35 j with EGS4 36 to model 
the electromagnetic interactions and PRESHOWER 
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code [37] for geomagnetic interactions. There is no 
significant dependence on the hadronic model because 
only photon-induced simulated showers are used in the 
method. The showers are simulated with thinning and 
the dethinning procedure is adopted [35] to simulate re- 
alistic shower fluctuations. 

The detector response is accounted for by using look- 
up tables generated by GEANT4 [39 simulations. Real- 
time array status and detector calibration information 
are used for each Monte Carlo (MC) simulated event. 
The Monte-Carlo events are produced in the same format 
as real events and analysis procedures are applied in the 
same way to both. The photon-induced MC set contains 
2 x 10 6 triggered events produced from 3380 CORSIKA 
showers by randomizing core location |40) . 



III. DATA SET 

We use the Telescope Array surface detector data set 
observed and recorded between 2008-05-f 1 and 2011-05- 
01. During this time period, the surface detector ar- 
ray was collecting data with a duty cycle greater than 
95% 0. 

We reconstruct each event with a joint fit of the geom- 
etry and Lateral Distribution Function (LDF) and de- 
termine the Linsley curvature parameter "a" (see Ap- 
pendix A for definition) along with the arrival direction, 
core location, and signal density at 800 meters S = Ssoo- 
As noted above, the same reconstruction procedure is 
applied to both data and Monte-Carlo events. 

For each real event "i" we estimate the energy of the 
hypothetical photon primary, El f — E J (S l ,8 l ,cf)' 1 ), i.e. 
the average energy of the primary photon, inducing the 
shower with the same arrival direction and S. The look- 
up table for E 7 (S,9,<fi) is built using the photon MC 
set; the dependence on azimuthal angle, tf>, is relevant for 
events with E 1 > 10 195 eV where geomagnetic preshow- 
ering is substantial. Photon-induced showers are natu- 
rally highly fluctuating. Consequently, the accuracy of 
the determination of E~ is about 50% at the one sigma 
level. In the present analysis, E 7 is used for event selec- 



tion only and therefore its fluctuations are well accounted 
for in the exposure calculation: the effect of these fluc- 
tuations is "lost" photons 30 , i.e. the photons with 
reconstructed energy below the energy cut. This will be 
estimated in Section Ivl 

We imposed the following requirements on both data 
and MC events: 

1. The shower core is inside the array boundary with 
the distance to the boundary larger than 1200 m; 

2. Zenith angle cut: 45° < 6 < 60°; 

3. The number of scintillator detectors triggered is 
>7; 

4. The joint fit quality cut, x 2 /d.o.f.< 5; 

5. S cut: E^(S l obs , 9\ <j) 1 ) > 10 19 eV or E y > 10 19 5 eV 
depending on the energy region discussed (the sec- 
ond variant is used for both Eq = 10 195 and 
E a = 10 20 eV). 

The cuts determine a photon detection efficiency which 
is greater than 50% for showers induced by primary pho- 
tons with energy above 10 19 eV. The calculation of expo- 
sure is given in Section [Vj The resulting data set contains 
877 events with E 1 > 10 19 eV and 45° < 9 < 60° which 
we used for our photon search. 

IV. METHOD 

To estimate the flux limit, we used an event-by-event 
method |30j . The Linsley curvature parameter "a" is 
used as a C-observable and S = S 800 is used as an E- 
observable. For each real event V we estimate the pair 
of parameters (S* bs , a l obs ) and the arrival direction (9 l , 
<fi l ) from the fit of shower front geometry and LDF. His- 
tograms of Linsley curvature are shown in Figure [2] 

We selected simulated gamma-induced showers com- 
patible with the observed 9 l , (j) 1 and Sl bs and calculate 
the curvature distribution of the simulated photon show- 
ers fl(a) as discussed in Reference [3U]. For each event, 
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FIG. 3: Left: Fit of the shower front for an event (2008-08-13 14:02:01, 6 =53.6°, E 1 = 1.29 x 10 19 eV, £=0.033) compared 
to average over all photon MC events with the same zenith angle and S. 68% of MC events have front delay within shaded la 
region. The front delay is counted from the plane front arrival time. Right: / 7 (a) for the same event; a f>s - observed value of 
curvature. The filled region indicates MC events with curvature smaller than a ^ s (3.3% of MC events). 



we determined the percentile rank of Linsley parameter 
a for photon primaries 



C = 



fUa)da, 



which is the value of the integral probability distribution 
function at the observed curvature. The shower front fit, 
/ 7 (a) and C for one of the events is shown in Figure [3] 

The distribution of C for data and MC is shown in 
Figure |ij Although the distribution of /* (a) varies with 
energy and arrival direction, C l for gamma-ray primaries 
would be distributed between and 1 uniformly by def- 
inition [J5]. O n the other hand, the actual distribution 
of C % in the data is strongly non-uniform (most of the 
events have C l below 0.5). 




FIG. 4: C distribution for the data set > 10 19 eV, 45° < 
9 < 60°. The black points show the data and the red line 
indicates the photon MC generated with and E~ 2 spectrum. 
The photon median is represented by the vertical gray line. 

Since the simulations of hadron-induced showers de- 
pend strongly on the hadronic interaction model, we do 



not use the hadronic showers simulations in calculation 
of the photon limit. 

Suppose that the integral flux of primary photons over 
a given energy range is F 1 . Then we expect to detect 



n(F 7 ) = (1 - X)F 7 A S 



(1) 



photon events on average, where A geom is the geometrical 
exposure of the experiment for a given data set and A is 
the fraction of "lost" photons (i.e. photons with primary 
energies within the interesting region which failed to en- 
ter the data set due to triggering efficiency and cuts). 

We calculate an upper limit on the primary photon flux 
based on the idea that photons satisfy a uniform distribu- 
tion from to 1 of the variable C. To do this, we examine 
all possible combinations of n events from the data set, 
where n covers the range from 3 to some large value M. 
We compare each combination to a uniform C distribu- 
tion using the Smirnov- Cramer- von Mises test [H], and 
let V(n) be the largest probability found in this way. By 
definition of the test V(0) = V(l) = V(2) = 1 and we 
assume M — 100 (for which all probabilities vanish in 
considered cases). See Appendix B for a description of 
the Smirnov-Cramer-von Mises test. To constrain the 
flux F~, at the confidence level £, we require 



M 



^P(n)W(n,n(F 7 ))<l-£, 



(2) 



where W(n,fi) is the Poisson probability of finding n 
events when the mean is n. To constrain the flux at 
95% confidence level (CL) we set £ = 0.95 and find h 
from The upper limit on the flux follows from ([I]), 



F 1 < 



(1 - \)A g 



(3) 



This method does not require any assumptions about 
hadron-induced showers and does not require the C- 
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57% 


70% 
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TABLE I: Relative exposure of TA SD (1 - A) to photons 
after consecutive application of cuts. 
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FIG. 5: Photon flux limits of the present work (TA) 
compared to the previous limits by AG AS A (A) [5], 
Yakutsk (Y) [7] and Pierre Auger Observatory (PA) [TT1 fl2] . 



observable to be strongly discriminating (like the muon 
density used in |7| [5] ) . 



V. EXPOSURE 

The geometrical exposure for the SD observation pe- 
riod with 45° < < 60° and boundary cut is 



igeojn 



1286 km sryr . 



(4) 



The fraction of the "lost" photons is calculated using 
a photon MC set generated with an E~ 2 spectrum. The 
values of (1 — A) after consecutive application of cuts are 
shown in Table |U 



VI. RESULTS 



Using the statistical method (Section IV I we arrive at 
the following results: 

n < 14.1 (95% CL), E y > 10 19 eV, 
n < 8.7 (95% CL), E 7 > 10 19 ' 5 eV , 
n < 8.7 (95% CL), E 7 > 10 20 eV. 



F 7 < 1.9 x 1CT 2 km" 2 sr- 1 yr" 1 (95% CL), E 1 > 10 19 eV, 
F 7 < 0.97 x 1(T 2 knr 2 sr _1 yr _1 (95% CL), E 7 > 10 195 eV, 
F 7 < 0.71 x 10~ 2 km" 2 sr- 1 yr _1 (95% CL), E 1 > 10 20 eV. 

These photon limits are shown along with the results 
of the other experiments in Figure [5] 

We obtain photon fraction limits by dividing the corre- 
sponding flux limits by the integral flux of the Telescope 
Array SD spectrum [18) : 

e 7 < 6.2% (95% CL), E 7 > 10 19 eV, 
e 7 < 28.5% (95% CL), £ 7 > 10 195 eV. 

The limits strongly constrain the top-down models of 
the origin of cosmic rays, see [32] for discussion. 

Next, we compare the results of event-by-event method 
with the results of simpler "photon median" method |11) . 
In the latter, the events having curvature greater than the 
median photon curvature are identified as photon candi- 
dates. This criteria corresponds to C > 0.5 and we ob- 
serve 3 candidate events with energy greater than 10 19 eV 
(see Figure [4} and no candidate events above 10 19 5 eV. 
This gives Poisson 95% confidence limits of n/2 < 8.25 
and n/2 < 3.09. The flux limits are F 7 < 2.3 x 10~ 2 , 
F 7 < 0.69 x 10~ 2 and < 0.51 x 10~ 2 km^sr-V -1 
for E a = 10 19 , 10 19 - 5 and 10 20 eV correspondingly. The 
limits using the two methods are in mutual agreement. 

Both the use of plastic scintillators sensitive to photon- 
induced showers and the application of event-by-event 
statistical method allowed us to put stringent limits on 
the flux of primary photons with energies in excess of 
10 eV with the data obtained during three years of the 
TA surface detector operation. The limits are strongest 
among those obtained in the northern hemisphere. The 
result depends neither on the choice of hadronic inter- 
action model, nor on possible systematics in the energy 
determination of hadronic primaries. 



Acknowledgments 

The Telescope Array experiment is supported by the 
Japan Society for the Promotion of Science through 
Grants-in-Aid for Scientific Research on Specially Pro- 
moted Research (21000002) "Extreme Phenomena in the 
Universe Explored by Highest Energy Cosmic Rays" , and 
the Inter-University Research Program of the Institute 
for Cosmic Ray Research; by the U.S. National Sci- 
ence Foundation awards PHY-0307098, PHY-0601915, 
PHY-0703893, PHY-0758342, and PHY-0848320 (Utah) 
and PHY-0649681 (Rutgers); by the National Research 
Foundation of Korea (2006-0050031, 2007-0056005, 2007- 
0093860, 2010-0011378, 2010-0028071, R32-10130); by 
the Russian Academy of Sciences, by the Russian Min- 
istry of Science and Education under contracts No. 8412, 
8525, by grants of the President of the Russian Federation 
NS-5590.2012.2, MK-1170.2013.2, Dynasty Foundation, 



6 



RFBR grants 11-02-01528, 13-02-01311 and 13-02-01293 
(INR), IISN project No. 4.4509.10 and Belgian Science 
Policy under IUAP VI/11 (ULB). The foundations of Dr. 
Ezekiel R. and Edna Wattis Dumke, Willard L. Eccles 
and the George S. and Dolores Dore Eccles all helped 
with generous donations. The State of Utah supported 
the project through its Economic Development Board, 
and the University of Utah through the Office of the Vice 
President for Research. The experimental site became 
available through the cooperation of the Utah School 
and Institutional Trust Lands Administration (SITLA), 
U.S. Bureau of Land Management and the U.S. Air Force. 
We also wish to thank the people and the officials of 
Millard County, Utah, for their steadfast and warm sup- 
port. We gratefully acknowledge the contributions from 
the technical staffs of our home institutions and the Uni- 
versity of Utah Center for High Performance Computing 
(CHPC). The cluster of the Theoretical Division of INR 
RAS was used for the numerical part of the work. 

Appendix A. LDF and shower front fit functions 

We perform joint fit of LDF and shower front with 7 
free parameters: x core , y core , 9, 4>, S 8 oo, h, a. 

S{r) = S soo x LDF(r) , 

*o(r) = to + t p i ane + ax 0.67 (1 + r I R L ) 15 LD F-°- 5 (r) , 

where t p i ane is a shower plane delay, a is a Linsley cur- 
vature parameter and the LDF(r) is defined as follows: 

LDF{r) =f(r)/f (800 m) , 

/ r v-i.3 / x -(.-1.2) ( 2X-0-6 



R m = 90 m, i?i = 1000 m, R L = 30 m, 
rj = 3.97- 1.79 x (sec(0) - 1) . 



Appendix B. Smirnov-Cramer-von Mises 
"omega-square" test implementation 

Let F(x) be theoretical distribution and F n (x) - ob- 
served distribution of n events. We define the distance 
between distributions by [4Tj : 

oo 

lo 2 = J (F n (C) - F{C)f dF(C) . 

— oo 

If C\ , C% , . . . , C n is a set of observed values in increasing 
order, lo 2 may be rewritten in the following form: 

2 1 " /2i-l ,„ A 2 

i=l v ' 
In this paper we compare distribution of event subset 

with uniform distribution U(0, 1). Therefore F(Ci) = Ci 

and we have: 




Required maximization of probability over subsets is 
therefore reduced to selection of n different events mini- 
mizing the above sum. The latter may be done with fast 
iterative procedure. 
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